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Abstract 

Second order perturbation theory and a Lipkin-Nogami scheme combined 
with an exact Monte Carlo projection after variation are applied to com- 
pute the ground-state energy of 6 < N < 210 electron-hole pairs confined 
in a parabolic two-dimensional quantum dot. The energy shows nice scaling 
properties as N or the confinement strength is varied. A crossover from the 
high-density electron-hole phase to the BCS excitonic phase is found at a 
density which is roughly four times the close-packing density of excitons. 
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I. INTRODUCTION 



As we understand, the interest in electron-hole states (excitons) in semiconductor physics 
is motivated by two facts. First, excitons have a bosonic character (as they are made up of 
a pair of fermions) and, thus, the many-exciton system is a candidate for a Bose condensate. 
This possibility was envisaged long ago |T[ , but regained attention in the last years after the 
Bose condensation of alcali atoms was achieved 0. Experimentally, signals of Bose-Einstein 
statistics have been identified in the photoluminiscence of quantum wells under strong laser 
pumping || , and indirect excitons in quantum wells are being manipulated via applied stress 
and inhomogeneous electric fields @j to reach the densities needed for Bose condensation. 
On the other hand, excitons are at the basis of many optical properties of semiconductors 
H . Recent experimental works have focused on the lowest dimensional structures, and very 
interesting properties have been found in the photoluminiscence of quantum wires || and 
single quantum dots [0. 

In the present paper, we study a two-dimensional quantum dot with a number of electron- 
hole pairs, 6 < iV < 210, i. e. intermediate between the very small dot [|7| and the quantum 
well 0. We study the dot at strong and intermediate confinement regimes by means of 
second-order perturbation theory and a variational (BCS) procedure. The main results of 
the paper may be summarised as follows. 

We found a breakdown of perturbation theory and a significant BCS pairing roughly 
at the same confinement strength, corresponding approximately to four times the close- 
packing density of excitons (or four times the Mott transition density). We notice that the 
BCS calculations were performed within the Lipkin-Nogami scheme |§ with exact projection 
onto the iV-pair sector || to avoid the incorrect behaviour of the naive BCS function in a 



finite system []10[ . Second, we found that the energy depends on N and the confinement 
strength in a scaled way. 



Our paper is complementary to |TT],[T2|, in which the multiexcitonic quantum dot is also 
studied. In [III], we found that the far- infrared absorption of the dot is dominated by a giant- 
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dipole resonance similar to the collective state apearing in nuclei |13| and metallic clusters 



|I4| . In paper ||12|| , the Bethe-Goldstone equations (the independent-pair approximation in 



Nuclear Physics [[15|]) are applied to study small (2 < iV < 6) clusters. The frequency for 
optical absorption with creation of an electron-hole pair shows a very interesting behaviour 
related to the apparent instability of the free (not confined) four-exciton cluster in two 
dimensions. 

II. PERTURBATION THEORY 

We study a direct-band-gap semiconductor with two parabolic bands. N electrons and 
iV holes are created by, e. g., strong laser pumping. We shall ignore recombination pro- 
cesses and electron-hole exchange. A model like this have been employed for the analysis 
of collective excitations in bulk semiconductors fllEf . The particles are forced to move in a 
two-dimensional region confined by a parabolic potential. This is a common approach in the 
study of self-assembled quantum dots RT7| . For simplicity, we take m h = m e and the same 
confining potential for both particles. Up to 210 electron-hole pairs will be allowed in the 
dot. 

The Hamiltonian of the system in oscillator units is written as 

TT 1 2N 

where uj is the dot frequency, q a = — 1 for electrons and +1 for holes. This Hamiltonian 



depends only on one constant, (3 = ^J(^^)/(huj) = yE c /(huj), where E c is the Coulomb 
characteristic energy, m is the electron effective mass, and k the dielectric constant of the 
material. In these units, the effective Bohr radius is a# = 1/(3. 

(3 —>■ is a high-density (strong confinement) limit in which the Bohr radius is much 
higher than the oscillator length (equal to one in our units). The independent-electron and 
-hole picture works in this limit, and the Coulomb interaction may be computed in pertur- 
bation theory. Notice that at high density, we have a system of independent fermions, not 
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bosons. This is one of the reasons preventing Bose condensation of excitons in homogeneous 
3D systems 

On the other hand, as (3 is increased, the dynamics become more and more dictated by 
the Coulomb forces. First, we shall observe the emergence of two-body correlations and 
the formation of electron-hole "Cooper" pairs (i. e. pairing in Fock space). With a further 
increase in (3, small excitons, biexcitons and higher complexes shall start playing a dominant 
role. 

Let us first consider the (3 — > limit, in which the Coulomb interaction may be computed 
in perturbation theory. We will make an additional simplifying assumption: N is such that 
there are N s heii closed shells in the (3 = limit, that is the number of electrons takes one of 
the following values N = N shell (N sheU + 1) = 6, 12, 20, 30, 42, . . . , 210. The ground state 
of such systems for small (3 values is spin-unpolarised, which means that both the electron 
and hole subsystems have total spin S = 0. The angular momentum of this state is L = 0. 

The (3 — > perturbative series take the form 

^ = b + b 1 (3 + b 2 p 2 + O((3 3 ), (2) 

where, the leading approximation to the energy is twice the energy of N independent elec- 
trons (or holes) 

b = 2NV4N + 1/3, (3) 
and for b\ and b 2 we arrive to the following expressions 

b\ = -2 (nijTixll/rlni.ni) -4 E ( n i> n 2 \l/r\n 2 , ni), (4) 

m<N/2 ni<n 2 <N/2 

b = _ 4 y- y- |£n<*/2( w 2,Ti|l/r-|n,ni)| 2 _ g (n 3 , n 3 \l/r\ ni , ni ) 2 

ni h/2n 2 h/2 e(n 2 )-e(n x ) ni h/2n :i h/2 <n 3 )-e(nl) 

- 24 T T (n3,n 3 \l/r\n 1 ,n 2 ) 2 _ 24 (n 3 , n 4 |l/r|ni, ni) 2 

ni <n 2 <N/2n 3 >N/2 2e (^) - e(nl) - e{n 2 ) ni < N /2 ni>n&>N , 2 <n 3 ) + e(n4) - 2e(n 1/ 

~ 8 E E {5(r2 3 ,n4|l/r|ni,n 2 ) 2 

ni<ri2<A r /2 n 4 >n 3 >A r /2 



+ 5(n 3 ,n 4 |l/r|n 2 ,ni} 2 - 4(n 3 , n 4 |l/r|ni, n 2 ) (n 3 , n 4 |l/r|n 2 , ni)} 

x (e(n 3 ) + e(n 4 ) - e(nl) - e(n 2 ))^ 1 . (5) 

The sums run over orbitals, which have been numbered sequentially. The sums over spin 
degrees of freedom have been explicitly evaluated. The Coulomb matrix elements are defined 

as 

(n 1 ,n 2 \l/r\n 3 ,n 4 ) = / -p ^^(^^(^^(rO^fa)- (6) 

The explicit form of the harmonic-oscillator orbitals is 

= ^l/l^'^r 2 )^ 272 ^, (7) 



where Ck,\i\ = Jk\/[n (k + |/|)!], and n = (k, I) is a composed index. The energy correspond- 
ing to <pk,i is e(k, I) — 1 + 2k + In terms of these energies, we have bo = 4Xi„ e i n )- 

Numerical values for the coefficients b\ and b 2 are presented in Table |. The sums entering 
the 62 coefficients were evaluated with a maximum of 20 shells. With respect to the number of 
shells included in the calculations, the convergence is slow, thus we used Shank extrapolants 
|18| to accelerate convergence. Fortunately, as a function of N, 6 2 saturates very fast and 
there is no need to perform calculations for N > 42. Notice the scaling laws 60 ~ §^V 3 / 2 , 
b\ ~ —0.96 iV 5 / 4 , 6 2 ~ —1-65 N for N > 42. b\ depends weaker on N (as compared with 
electrons, for which the power is 7/4 instead of 5/4 Jl9|) because of the partial cancellation 
between attractive and repulsive Coulomb matrix elements. 

Least-squares fits to the data in Table | lead to 

* -0.960853- + (8) 



AT5/4 - N N 2 N 3 

iV = - L665 + ^\T-^- (9) 



A. Approximate scaling and breakdown of the perturbative expansion 



For large enough N, we can use the asymptotic expressions for the coefficients to show 
that E/N 3 / 2 is approximately a function of the combination iV _1//4 /3. 
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The physics behind the scaling (10) is the following. As a result of cancellation between 
Coulomb attraction and repulsion, the size of the system is practically constant. Then, by 
increasing N, we increase the density and depress the effects of the Coulomb interaction. 
Approximate scaling of the energy is also characteristic of confined electron systems |H| and 
charged bosons in two dimensions [f20fl . Notice that in a pure electron system, where the 
interparticle potential is always repulsive, an increase in N leads to a decrease of the density 
and an enhancement of correlation effects. The energy turns out to be a function of iV 1//4 /? 
at low (3. 

We shall stress that the scaling law ([H]) is expected to be observed also in the strong 
coupling, (3 — > oo, limit in which the energy shall be roughly proportional to the energy of 
N independent excitons, 



E 

Tiuj 



= a (3 2 + ---, (11) 



where a ~ —N. The right hand side of Eq. ([11]) may thus be written as N 3 ^ 2 (—P 2 /N 1 ^ 2 ). 
The variational results of the next sections also support the scaling behaviour (|10"1). 

A naive estimation of the convergence radius for the series (0) gives j3 < (3 C = bi/b 2 . This 
estimation may be obtained formally as the pole of the Pade approximant 

PxM = bo + -r^-s, qi = -b 2 /h, (12) 
1 + qiP 

which reproduces the expansion (0) for small (3 values. Notice the high-iV asymptotic be- 
haviour, (3 C ~ 0.58 iV 1//4 . (3 C gives an estimate for the density at which exciton effects 
become important. Indeed, the density in our units is p ~ N/(ir(r 2 )) w 3 iV 1/,2 /(27r), thus 
f3 c may be expressed in terms of p. Turning back to ordinary units, we get a critical density, 
p c ~ 1/(0.24 Ti a 2 B ), i. e. approximately four times the close-packing density of excitons, 
l/(7r a 2 B ). This fact is consistent with the belief that screening is less effective in two dimen- 
sions. Below p c , exciton effects shall dominate the quantum dynamics. As will be seen, p c 
is also at the onset of pairing in the BCS estimate of the next sections. 



In the following sections, we will perform variational estimations expected to be valid 
when pairing is not so strong, that is in the regime /3/iV 1//4 < 1. 

III. VARIATIONAL ESTIMATIONS 

Let us turn to the variational calculations. The simplest variational estimation one can 
try is first-order perturbation theory. 

^<E PT1 ((3) = b + b 1 (3. (13) 

This estimate may be improved by introducing a frequency, Q, as an additional varia- 
tional parameter, i.e. by taking as trial function the product of two Slater determinants of 
harmonic-oscillator states with a frequency Q. The result is, 

A <minn {l (f! + 1/ 0) £pT1 ( s ^ /3 )}. (14) 

We checked that the result coming from ( |T4"D practically coincides with the Hartree-Fock 



(HF) energy for this system [nj. Thus, we will call flTjj ) the HF estimate. 

The mechanism by which the energy is lowered is pairing. We may take account of it 
with the help of a BCS-like wave function This may be a good estimation for weak 



pairing, when correlations are not so strong. In the (3 axis, it means < 1. The wave 

function is given by 

max 

\BCS) = ]J K' + Vjh+e}) \0) h |0) e . (15) 
i=i 

hj and elj are hole and electron (harmonic oscillator) creation operators acting on their 
respective vacua 10)^ and |0) e . j = (k,l,s z ) is a composed index, f = (k, —I, —s z ). s z is 
the spin projection. Vj and Uj are normalised according to u^j + 1>| = 1. The total angular 
momentum corresponding to \BCS) is zero because the angular momentum of each pair is 
zero. The mean value of the total electron (hole) spin may be forced to be zero by requiring 
v(k, I, s z ) = v(k, I, —s z ). Thus, Vj does not depend on s z and we can write v n instead of Vj. 
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\BCS) is not an eigenfunction of the particle number operator. In a finite system, we 
shall project onto the state with the correct number of particles. This will be done in two 
steps: first, an approximate projection before variation over the parameters v n entering the 
BCS function (the Lipkin-Nogami scheme ||), and then an exact Monte Carlo projection 
of the BCS function onto the sector with N pairs ||. 



A. The Lipkin-Nogami estimate 

In the Lipkin-Nogami (LN) method ||, one assumes an approximate polynomial depen- 
dence of H on the particle number operator N, 

H = X + 2X 1 N + X 2 N 2 . (16) 

By taking expectation values of H over exact and BCS functions and comparing results, 
we arrive to 

E LN = E BCS - 2Ax ((N)bcs -N)-X 2 ((N 2 ) bcs - N 2 ) , (17) 

where 

Ebcs = (H)bcs = E{ 4e " - 2P(n,n\l/r\n,n)}v 2 

n 

-2(3 (ni,n 2 \l/r\n 2 ,n 1 ){v 2 1 v 2 2 +v ni u ni v n2 u n2 }. (18) 



Minimization over Aj leads to 



N = (E e^acs = (E hjh^Bcs 
j j 

= 2E <■ (19) 



The equation of minimum with respect to the v n can be written in the form of standard 
jap equations 

A n = /3 y (n,7ii|l/r|ni,n) . Ani (20) 
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where the HF energies are given by 



et? F = e n — — (n, n\l/r\n, n) - (3 ^ (n, ni|l/r|ni, n) v\ x — X 2 (N — v 2 ), (21) 



2 



and we used the common BCS parametrization 



1 ■ , (22) 

The chemical potential, /x = A1 + A2/2 was introduced in Eq. (p^). For the determination 
of A 2 , the system of equations 

(H - A - 2AxiV - A 2 2V 2 ) BC5 = 0, (23) 
((H-X -2X 1 N-X 2 N 2 )N) BC s = 0, (24) 
((# - A - 2AxiV - X 2 N 2 )N 2 ) BCS = 0, (25) 

is used |5| . It makes the LN method not throughly variational. The first equation determines 
the constant Aq. The second turns to be equivalent to the gap equation (GO). For A 2 , we get 



aid?, — a 2 a^ 

X 2 = — , (26) 

^305 — a 2 

where 

Oi = (HN 2 ) BCS - (H) BCS (N 2 ) BCS , (27) 
a 2 = (N 3 ) BCS - N(N 2 ) BCS (28) 
a 3 = (HN*) BCS - (N 2 ) 2 BCS1 (29) 
a 4 = (^iV) B cs - N{H) BCS , (30) 

« 5 = (iV 2 ) B cs - N 2 . (31) 

The resulting equations were solved iteratively starting from e^ F = e n , A n = 0.2. First, 
the explicit expressions for v 2 are used and the nonlinear equation ([19]) is solved for /j. After 
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that, we obtain A2 from (p6[), and the A n and e^ F are recalculated from ( |20| , [2~T1) . The process 
is repeated until the variation in any of the e^ F is less than 1CT 10 . 

Calculations were carried out for 6 < N < 90 pairs and a maximum of 600 one-particle 
states for both electrons and holes (i. e. 300 orbitals, because there are 2 spin states for 
each orbital). The absolute error in computing Coulomb matrix elements is less than 10~ 8 . 
As is Eq. (|i~4l) , we introduced an additional parameter Q, and used the inequality 

E< mi nJ±(n + l/n)E L j£^p\\, (32) 



where Eln is the result from Eq. (|17|) at Q — 1. The variation of f2 can be thought of as 
a simplified self-consistent Hartree-Fock-Bogoliubov procedure, in which the mean field is 
forced to be a harmonic potential. 

We show in Fig. p] the energy coming from Eq. (^) versus j3 for N = 6 and 90 pairs (the 
curves LN). The HF estimates (|14]) and BCS curves (the A2 = limit of LN) are also given 
for comparison. The lowest curves, labelled "proj", correspond to the exact projection of the 
next section. We notice that the LN method moves the BCS threshold for pairing towards 
zero, but the energy itself remains very close to the HF curve. A significant departure occurs 
only for (3 > (3 C . {(3 C « 0.86 for iV = 2, and 1.7 for N = 90). 

B. The Monte Carlo projection 

The next step is an exact projection of the wave function onto the iV-pair sector. The 



situation is similar to the calculations carried out for nuclei |p2jl , where a variety of projection 
methods have been developed. 

We project the wave function after the {v n } are determined for given f3 and Q. With 
this function, the mean value of the hamiltonian at a shifted (3 is computed and multiplied 



by the factor given in Eq. fl32|). The projected energy takes the following expression, 

E \n= E in./'. ,/vS (./, ,/v). (33) 
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where the j's have the same meaning as in Eq. Ql5|), and the sum runs over possible 
combinations of N states, {ji, . . . ,Jn}, from a maximum of N max states allowed in the LN 
calculation. The "weights", W, and "energies", e, are defined as 

v 2 v 2 ( v 2 ...vl V 1 

(./: ./v) 2 2 ^--z 3 E dJ\i/r\jj)-P E 0',/|iAb",j> 

j'e{ii,...,ijv} ie{ji,...,ijv} j',i'e{ii,-,ijv}> jW 

E (./-./' l/r./'../}^ 1 . (35) 

The expression (^) for the projected energy allows a simple Monte Carlo evaluation, 
where the sets . . . ,j^} are generated with probability W(ji, . . . , jx) by means of a 
Metropolis algorithm [KB]. Other equivalent forms of Eq. Q3"5|), see for example Ref. P3 



arc 



not suited for this evaluation. The procedure seems to be particularly efficient in Nuclear 
Physics calculations as well Q. 

The results are also drawn in Fig. |l|. The improvement is significant for (3 ~ /3 C , and its 
relative importance diminishes as N is increased. 

C. Results 

We show in Fig. |^ our best results for the energies of the systems under study. The 
scaled energies show a remarkable similarity. A significant pairing (i. e. departure from the 
HF curve) is seen only for > 0.55. 

Finally, we give a parametrisation of the ground-state energy obtained from the best of 
our variational estimates. The energy is written in the form of a Pade approximant |2B , 



E *- = b ° + W+ l + qiP + q2 f? > ( 36 ) 

where p% = qlpA/q2 — blq2, and the coefficients p±, qi, and q2 are fitted from our numerical 
results. The obtained values are shown in Table |. 
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IV. CONCLUDING REMARKS 



We have studied electron-hole systems in quantum dots under strong and intermediate 
confinement, where the dense electron-hole or the BCS excitonic phases are present. 

The breakdown of perturbation theory and a significant pairing in the BCS wave function, 
both take place at a density which is roughly four times the close-packing density of excitons. 
We interpret this result as a crossover between the two phases. 

As mentioned before, with an increase in f5, particle correlations shall become more and 
more important. We shall observe signals of the "excitonic", "biexcitonic" , etc. insulat- 
ing phases. The true energies and wave functions of these phases shall be obtained with 



more powerful methods as, for example, Green-function Monte Carlo method p3fl . Even 
a variational Monte Carlo estimation, as that one carried out for the homogeneous case in 



Ref. p6| , may be biased by the chosen trial functions. The density matrix renormalisation 
group method of Ref. ]27] could also be useful. The very interesting question about whether 
the system remains bound after the external potential is switched off, still remains to be 
answered. We have some indications that the two-dimensional triexciton is bound and the 



four-exciton system is unbound fl2 |. But the situation may be analogous to nuclei, where 
there is a small instability island around atomic number 5. Some of these problems are 
currently under investigation. 
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FIGURES 

FIG. 1. a) and b): Ground-state energies of the 6-exciton and 90-exciton systems respectively. 



FIG. 2. Scaling of the ground-state energies 
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TABLES 

TABLE I. The coefficients b\, 62, P3, P4, Qi, and q 2 . 



N 


61 /N 5 / 4 


b 2 /N 


P3 


Pa 


<7i 


92 


6 


-1.001 


-1.58 


21.9817 


-461906 


15.3098 


9.92132 


12 


-0.983778 


-1.62 


43.4125 


-47.3284 


6.46992 


4.84953 


20 


-0.9758 


-1.64 


53.9705 


-47.5022 


3.60744 


2.85761 


30 


-0.971391 


-1.65 


84.1615 


-64.3801 


3.68693 


2.58195 


42 


-0.968681 


-1.65 


127.946 


-115.014 


6.85897 


3.44577 


56 


-0.966881 




172.825 


-157.173 


8.02801 


3.55984 


72 


-0.965621 




213.377 


-189.426 


8.25862 


3.35560 


90 


-0.9647 




259.185 


-224.076 


8.27331 


3.16176 


110 


-0.964009 












132 


-0.963475 












156 


-0.963055 












182 


-0.962712 












210 


-0.962433 
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